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ABSTRACT:  Dislocation  mediated  plasticity  in  the  a  phase  of  the  energetic 
molecular  crystal  cyclotrimethylene  trinitramine  (RDX)  was  investigated  using  a 
combination  of  atomistic  simulations  and  the  Peierls— Nabarro  (PN)  model.  A 
detailed  investigation  of  core  structures  and  dislocation  Peierls  stress  was 
conducted  using  athermal  atomistic  simulations  at  atmospheric  pressure  to 
determine  the  active  slip  systems.  Generalized  stacking  fault  energy  surfaces 
calculated  using  atomistic  simulations  were  used  in  the  PN  model  to  also 
estimate  the  critical  shear  stress  for  dislocation  motion.  The  primary  slip  plane  is 
found  to  be  (010)  in  agreement  with  experimental  observations,  with  the  (010) 

[100]  slip  systems  having  the  lowest  Peierls  stress.  In  addition,  atomistic 
simulations  predict  the  (02l)[0l2],  (02l)[l00],  (01l)[l00],  (00l)[l00],  and 
(00l)[010]  slip  systems  to  have  Peierls  stress  values  small  enough  to  allow 
plastic  activity.  However,  there  are  less  than  five  independent  slip  systems  in  this  material  in  all  situations.  The  ranking  of  slip 
systems  based  on  the  Peierls  stress  values  is  provided,  and  implications  are  discussed  in  relation  to  experimental  data  from 
nanoindentation  and  shock-induced  plastic  deformation. 


1.  INTRODUCTION 

Cyclotrimethylene  trinitramine  (RDX),  an  energetic  molecular 
crystal,  is  used  in  many  explosive  formulations  for  both  civilian 
and  military  purposes.  Initiation  and  detonation  in  explosives  is 
considered  to  happen  by  formation  of  regions  of  localized 
energy  called  “hot  spots.”  Theoretical  considerations1  have 
proposed  that  hot  spots  need  to  have  micrometer  dimensions 
and  sustain  high  temperatures  (~500  °C)  for  microseconds. 
Various  theories  have  been  proposed  regarding  the  precise 
mechanisms  of  hot-spot  formation.  In  polymer  bonded 
explosives,  local  heterogeneity,  presence  of  defects  (disloca¬ 
tions,  pores,  cracks),  and/or  adiabatic  heating  due  to  plastic 
deformation  are  thought  to  contribute  toward  formation  of  hot 
spots.  Experimental  and  theoretical  studies  in  single  crystals 
have  suggested  different  mechanisms.  Dislocation  pileup  and 
subsequent  energy  release  on  depinning  have  been  proposed2 
as  a  mechanism  for  hot-spot  formation  in  single  crystals. 
Anisotropic  activation  of  slip  systems  was  found  to  correlate 
with  the  observed  anisotropic  response  of  shock-loaded 
pentaerythritol  tetranitrate  single  crystals,  and  this  led  to  the 
development  of  the  steric  hindrance  model.3  According  to  this 
model,  the  lack  of  activated  slip  systems  for  a  particular 
crystallographic  orientation  leads  to  deformation  of  molecules, 
resulting  in  bond  breaking  and  subsequent  release  of  energy 
which  eventually  produces  initiation. 

Owing  to  their  suspected  role  in  processes  leading  to 
detonation,  plastic  deformation  mechanisms  operating  in  RDX 
have  been  investigated  by  computational  and  experimental 

ACS  Publications  ©2013  American  Chemical  Society  5326 


means.  Dynamic  yield  points  determined  from  shock  loading  of 
single  crystals  at  a  Rankine— Hugoniot  pressure  (Rrh)  of  2.25 
GPa  for  (ill),  (021),  and  (100)  planes  were  reported  to  be 
isotropic.4  The  response,  described  in  terms  of  the  shape  of  the 
velocity  history  profiles,  was  found  to  exhibit  significant 
variation  depending  on  the  shock  direction.  Plate  impact 
experiments  in  the  (ill)  direction5  indicated  a  change  in 
response  as  the  pressure  varied  from  1.25  to  2.25  GPa.  This 
variation  was  attributed  to  the  nucleation  of  (00l)[010]  partials 
upon  shock  loading  and  subsequent  plastic  hardening,  which 
was  validated  using  molecular  dynamics  simulations.  The 
observation  was  also  used  to  explain  the  anomalous  response 
of  (021)  shocked  RDX  crystals.6  High  pressure  (7—20  GPa) 
plate  impact  experiments  on  the  (ill),  (210),  and  (100) 
planes7  revealed  considerable  anisotropy  in  shock  decom¬ 
position  of  RDX.  Crystals  shocked  perpendicular  to  the  (210) 
and  (lOO)  planes  were  found  to  be  more  sensitive  compared  to 
those  shocked  perpendicular  to  the  (ill)  plane.  In  the  present 
study  we  propose  an  interpretation  of  this  observation  based  on 
the  mechanisms  of  dislocation-mediated  plasticity. 

X-ray  diffraction,  microindentation,  and  etch-pitting  were 
used  to  correlate  deformation  features  with  possible  active  slip 
systems  in  RDX.  Experimental  studies  from  different  sources 
have  proposed  (010)  as  the  primary  slip  plane.8,9  Micro- 
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Figure  1.  (a)  Unit  cell  of  cr-RDX  in  the  [100]  projection,  (b)  a-RDX  crystal  in  the  (100)  projection  with  the  slip  planes  considered  in  this  study 
shown.  The  two  types  of  slip  plane  (PI  and  P2)  are  indicated  for  (001). 


indentation  studies  by  Gallagher  et  al.10  suggested  slip  systems 
of  the  {021 } [  100]  and  (010)[00l]  type.  Investigations  by 
Halfpenny  et  al.11  also  suggest  (010)  as  the  primary  slip  plane 
and  hinted  at  the  possibility  that  {Oil}  or  {021}  could  be 
secondary  slip  planes.  Nanoindentation  experiments  performed 
on  (210),  (021),  and  (001)  oriented  crystals  of  RDX  showed 
elastic— plastic  load  excursion  at  estimated  resolved  shear 
stresses  within  the  G/15  to  G/10  range,12  where  G  is  the 
shear  modulus.  These  observations  were  attributed  to  the 
homogeneous  nucleation  of  dislocations.  The  hardness  values 
were  reported  to  be  nearly  isotropic,  but  significant  variation 
was  observed  for  the  load  at  the  first  “pop-in”  event. 
Tomography  of  indent  impressions  revealed  features  consistent 
with  slip  on  (010),  {021},  and  {011}  planes,  and  it  was 
conjectured  that  RDX  exhibits  hindered  plasticity  due  to  the 
lack  of  independent  slip  systems.  Munday  et  al.13  investigated 
the  competition  between  ductile  deformation  and  brittle 
fracture  using  generalized  stacking  fault  surfaces  coupled  with 
Rice’s  criterion  for  dislocation  emission  from  the  crack  tip.14 
They  found  that  the  (010),  (Oil),  (021),  and  (001)  planes 
possess  slip  systems  exhibiting  ductile  behavior. 

A  consensus  exists  in  the  literature  regarding  slip  on  the 
(010),  {021},  and  {011}  planes.8'9'11-13  Recent  molecular 
dynamics  simulations  have  also  indicated  possible  slip  on  (001) 
planes.5'13  Activation  of  slip  systems  requires  that  the  resolved 
shear  stress  is  higher  than  the  critical  stress  (Peierls  stress) 
required  for  the  motion  of  individual  dislocations  on  the 
respective  slip  system.  Knowledge  of  Peierls  stresses  is  essential 
in  understanding  preferential  activation  of  slip  systems  under  a 
general  state  of  stress. 

In  this  study,  fully  atomistic  simulations  are  performed  to 
evaluate  this  critical  stress.  The  investigation  was  conducted  for 
the  (010),  (001),  (Oil),  and  (021)  crystallographic  planes  in  cr- 
RDX,  which  are  shown  in  Figure  1.  As  discussed  above,  these 
planes  are  expected  to  display  ductile  behavior.  The  [  100]  slip 
direction  in  all  the  above  planes  is  studied,  which  has  been 
proposed12’13'15  as  a  possible  cross-slip  direction.  In  addition, 
the  (010) [001],  (001) [010],  (02l)[0l2],  and  (01l)[0ll] 
systems  are  considered,  which  have  been  also  suggested  in 
the  literature  as  possible  slip  systems.13,15 


2.  METHODOLOGY 

cr-RDX  is  modeled  using  the  flexible  Smith— Bharadwaj 
potential  developed  for  condensed  phase  nitramines.16  The 
potential  has  been  reported  to  accurately  reproduce  elastic  and 
thermal  properties  of  RDX.5  Nonbonded  interactions,  which 
includes  electrostatic  and  van  der  Waals  interactions,  are 
explicitly  considered  up  to  a  cutoff  radius  of  10  A. The  charges 
are  modified  to  account  for  the  effects  of  polarization  in  the 
condensed  phases.16  Bonded  interactions  are  represented  with 
angle-bending,  bond-stretching,  and  dihedral  interactions. 
Compared  with  the  other  potential  available  for  RDX,  the 
Sorescu— Rice— Thompson  potential,17  which  represents  the 
molecules  as  rigid  entities,  the  Smith— Bharadwaj  potential 
allows  molecular  flexibility.  As  discussed  in  this  paper, 
molecular  conformation  changes  are  important  in  the 
mechanics  of  crystal  defects.  This  has  determined  the  choice 
of  force  field  in  this  work.  A  schematic  of  the  simulation  setup 
is  shown  in  Figure  2.  For  the  (010)  and  (001)  planes,  the 
simulation  cell  is  of  size  100  X  2  X  100  unit  cells  respectively  in 
the  n,  p,  and  q  directions  and  consists  of  3  360  000  atoms.  The 
simulation  cell  size  was  selected  as  large  as  possible  in  order  to 


Figure  2.  Simulation  setup.  In  the  figure,  p,  q,  and  n  form  an 
orthogonal  coordinate  system.  In  all  simulations,  the  glide  plane 
normal  (g)  is  aligned  with  n.  For  edge  dislocation  dipoles,  the 
dislocation  line  direction  (|)  is  aligned  with  p  and  the  Burgers  vector 
(b)  is  aligned  with  q.  For  screw  dislocation  dipoles,  |  and  b  are  aligned 
along  p.  ra  represents  the  applied  shear  stress  in  the  b  direction  and  is 
shown  in  the  figure  as  applied  to  an  edge  dislocation  dipole.  L  denotes 
the  length  of  the  simulation  cell. 
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reduce  the  effect  of  the  boundary  conditions  on  the  quantities 
of  interest,  while  staying  within  our  current  computational 
capabilities.  In  all  simulations,  the  glide  plane  normal  (g)  is 
aligned  with  n.  For  edge  dislocation  dipoles,  the  dislocation  line 
direction  (§)  is  aligned  with  p  and  the  Burgers  vector  (b)  is 
aligned  with  q.  For  screw  dislocation  dipoles,  £  and  b  are 
aligned  with  p.  Long-range  interactions  are  calculated  using  the 
particle-particle-particle-mesh  (PPPM)  algorithm. 

For  the  (Oil)  and  (021)  planes,  periodicity  in  the  glide  plane 
requires  use  of  nonorthogonal  simulation  cells.  For  simulations 
on  these  slip  planes,  the  simulation  cell  were  of  sizes  30  X  2  X 
50  and  20  X  2  X  32,  respectively,  in  the  n,  p,  and  q  directions. 
The  Ewald/N  algorithm  was  used  to  calculate  long-range 
interactions  in  simulations  with  nonorthogonal  cells.  In  all 
cases,  a  vacuum  pad  is  included  in  the  direction  perpendicular 
to  the  glide  plane.  Both  orthogonal  and  nonorthogonal  cells  are 
shown  in  Figure  3. 


Figure  3.  Schematic  of  the  simulation  cells  used  in  the  study.  Z 
represents  the  simulation  box  vector  used  to  define  the  periodicity,  n, 
q,  and  L  are  defined  in  Figure  2.  (a)  Orthogonal  cells  were  used  for  the 
(010)  and  (001)  slip  planes.  (b)  Nonorthogonal  cells  were  used  for 
the  (Oil)  and  (021)  slip  planes. 


a-RDX  has  an  orthorhombic  structure  with  the  space  group 
Pbca.  There  are  eight  molecules  per  unit  cell,  with  all  molecules 
in  the  Caae  conformation.  The  conformational  state  of  RDX 
molecule  is  defined  by  the  shape  of  the  ring  (chair  denoted  by 
“C”,  boat  denoted  by  “B”,  or  twist  denoted  by  “T”)  and  the 
position  of  the  three  nitro  groups  (axial  denoted  by  “a”  and 
equatorial  denoted  by  “e”).  A  nitro  group  oriented  in  the 
direction  parallel  to  the  ring  normal  is  considered  axial. 

The  initial  structure  was  created  using  neutron-diffraction 
data  from  Choi  and  Prince19  and  was  allowed  to  relax  at  0  K 
and  atmospheric  pressure.  A  dislocation  dipole  was  created  by 
imposing  the  Volterra  displacement  field  on  the  centers  of  mass 
(COM)  of  the  molecules,  which  ensures  that  the  molecular 
bonds  and  angles  are  not  distorted  by  the  applied  displacement 
field.  The  individual  dislocations  of  the  dipole  are  arranged  on 
the  same  glide  plane  in  such  a  way  that  a  periodic  array  of 
dipoles  results  (on  the  glide  plane)  upon  application  of  periodic 
boundary  conditions.  The  distances  between  the  dislocations  in 
the  unit  cell  and  that  between  these  and  their  images  are  equal. 
This  ensures  that  the  Peach— Koehler  force  acting  on  each 
dislocation  due  to  the  interaction  with  others  in  the  periodic 
array  vanishes  in  the  initial,  unloaded  configuration.  This  holds 
in  materials  with  generic  anisotropy  since  the  resolved  shear 
stress  in  the  direction  of  the  Burgers  vector  and  in  the  glide 
plane  is  always  symmetric  relative  to  the  position  of  the  core.18 
An  image  exists  in  the  direction  of  the  vacuum  pads,  but  no 
infinite  array  of  images  exists  in  the  n  direction.  When  the 
model  is  distorted  and  the  dislocation  cores  shift,  the 


interaction  with  images  in  the  glide  plane  produces  a  net 
Peach— Koehler  force  on  each  dislocation  of  the  dipole.  This 
introduces  a  small  error  in  the  estimation  of  the  Peierls  stress. 
This  error  was  estimated  using  the  linear  elastic  solution,  based 
on  the  shift  of  the  cores  by  b/2,  i.e.,  a  shift  corresponding  to  the 
peak  of  the  Peierls  barrier,  and  was  found  to  be  less  than  1%  in 
all  cases  considered. 

It  should  also  be  noted  that,  since  each  molecule  remains 
electrically  neutral  during  deformation,  the  cores  do  not  carry 
electric  charges.  Electric  dipoles  may  exist  due  to  the  large 
deformations,  but  their  field  decays  faster  than  the  elastic  field; 
hence  the  dipole— dipole  interaction  of  the  cores  with  their 
images  is  weak. 

An  infinitely  long  (and  straight)  dislocation  is  created  upon 
application  of  periodic  boundary  conditions  in  the  p  direction. 
A  resolved  shear  stress  (r)  is  imposed  so  that  the  individual 
dislocations  of  each  dipole  move  toward  each  other,  i.e.,  toward 
the  center  of  the  simulation  cell.  This  is  done  by  applying  a 
homogeneous  shear  strain  to  all  COM  of  molecules  in  the 
relaxed  model,  followed  by  relaxation  while  holding  the 
boundaries  of  the  model  fixed.  The  boundary  is  maintained 
fixed  by  zeroing  the  forces  in  all  directions  on  a  layer  of  two 
unit  cell  thickness  at  the  boundary.  The  presence  of  the  vacuum 
pads  in  the  n  direction  facilitates  the  application  of  this 
distortion.  To  prevent  significant  steric  hindrance  due  to  the 
rigid  shift  of  the  molecules  in  this  step,  the  relaxation  is  divided 
in  two  stages.  In  the  first  stage  a  repulsive  potential  is  added  to 
the  nonbonded  interactions  to  help  separate  the  potentially 
overlapping  molecules  and  few  relaxation  steps  are  taken.  Then, 
this  dummy  field  is  removed  and  relaxation  proceeds  normally. 
In  addition,  the  displacement  imposed  in  each  loading  step  is 
kept  small  to  reduce  the  likelihood  of  overlaps. 

The  applied  shear  strain  was  increased,  starting  each  time 
with  the  relaxed  stress-free  dislocation  structure,  until  the 
dislocation  was  found  to  shift  by  at  least  one  lattice  spacing. 
The  corresponding  stress  value  is  the  Peierls— Nabarro  stress 
(rPN)  of  an  infinitely  long,  straight  dislocation  of  the  particular 
type.  For  partials,  the  Peierls— Nabarro  stress  is  calculated  as  the 
sum  of  the  critical  stress  and  the  contribution  to  the  Peach— 
Koehler  force  from  the  stacking  fault.  When  considering 
partials,  the  stacking  fault  is  always  located  between  the  two 
dislocations  explicitly  represented  in  the  model;  hence  its 
action  is  to  pull  them  toward  each  other.  The  stacking  fault 
energy  was  obtained  from  the  literature.13 

All  simulations  were  performed  using  LAMMPS,  a  molecular 
simulation  code  from  Sandia  National  Laboratories.20  Mini¬ 
mizations  were  performed  using  the  conjugate  gradient 
algorithm.  Residual  kinetic  energy  was  drained  by  running 
damped  dynamics  with  a  damping  constant  of  20  kcal/mol-fs 
until  the  temperature  of  the  simulation  cell  was  below  0.2  K. 
Simulations  were  performed  using  1024  processors  on  the 
Blue-Gene/L  at  the  Computational  Center  for  Nanotechnology 
Innovations  (CCNI)  at  Rensselaer  Polytechnic  Institute. 

3.  RESULTS 

3.1.  Structure  of  Dislocation  Cores.  The  core  structure  is 
evaluated  by  calculating  the  relative  displacement  of  the  COM 
of  molecules  across  the  glide  planes  («(»)). The  width  of  the 
dislocation,  w,  is  defined  as  the  distance  between  the  two 
inflection  points  of  u(x),  i.e.,  approximately  b/ 4  <  u(x )  <  3F/4. 
A  wide,  planar  core  would  be  easier  to  move  than  a  compact 
core.  The  equation  u(x )  =  b  arctan[(x  —  x0)/iv]  is  fitted  to  the 
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Figure  4.  Structure  of  the  (010)  [100]  edge  dislocation,  (a)  The  dislocation  with  b  =  [100]  splits  into  two  partials  of  b  =  l/2[  100] .  (b)  Structure  of 
the  b  =  1/2[100]  partial. 


core  profiles  of  stable  dislocations  resulting  from  the  simulation 
in  order  to  evaluate  the  widths. 

Dislocations  are  considered  to  be  stable  if,  upon 
minimization,  the  molecular  layers  on  either  side  of  the  glide 
plane  did  not  separate  (cleave)  and  the  magnitude  of  the 
Burgers  vector  remained  equal  to  that  corresponding  to  the 
applied  displacement  field.  The  various  slip  planes  considered 
in  the  study  are  schematically  represented  in  Figure  1,  and  the 
core  structures  obtained  from  atomistic  calculations  are 
presented  in  this  section.  For  each  crystallographic  plane 
there  exist  two  parallel  slip  planes  which  differ  in  the  nature  of 
intermole cular  bonding  across  the  glide  plane.  An  example  is 
shown  in  Figure  1,  where  it  is  apparent  that  two  different  types 
of  (001)  planes  form  the  stacking,  with  two  others  in  the 
sequence  being  identical,  but  rotated  by  180°  relative  to  the 
[001]  axis.  The  planes  differ  in  attachment  (cohesive)  energies 
as  described  previously.13  In  the  present  study  only  the  planes 
with  lower  cohesive  energies  are  considered. 

3.1.1.  (010)  Plane.  Two  directions  of  the  Burgers  vector  are 
considered  in  the  (010)  plane:  the  [100]  and  the  [001] .  In  each 
case  the  screw  and  edge  dislocations  are  studied.  The  (010) 
[100]  slip  system  has  been  reported  in  previous  studies  to 
exhibit  dislocation  activity.12,13  Slip  was  also  reported  in  the 
(010)[00l]  system  in  previous  experimental  investigations.8,10 
The  [100]  direction  was  reported13  to  have  a  stable  stacking 
fault  at  b  =  1/2[100]  with  a  stacking  fault  energy  of  101  mj/m2. 

Let  us  consider  first  the  edge  dislocation.  It  was  found  in  our 
simulations  that  the  edge  dislocations  of  the  (010)[l00]  type 
with  b  =  [100]  are  not  stable  and  dissociate  into  partials  of  pure 
edge  character  with  b  =  1/2[100],  therefore  providing  support 
to  the  observation  of  a  minimum  in  the  [100]  trace  of  the  y- 
surface  at  a  relative  shift  of  1/2  [  100] .  The  structure  of  the 
relaxed  core  is  compared  to  the  structure  of  the  Volterra 
dislocation21  in  Figure  4a.  Figure  4b  shows  the  structure  of  a  b 
=  1/2[100]  partial  dislocation.  The  split  of  the  dislocation  with 
b  =  [100]  in  two  partials  is  clearly  observed  in  Figure  4a,  while 
the  partial  dislocation  with  b  =  1/2[100]  has  a  compact  core, 
very  close  to  that  of  the  Volterra  dislocation. 

Full  screw  dislocations  of  the  (010)[l00]  type  were  stable 
when  t  =  0,  but  were  found  to  split  into  partials  with  b  =  1/ 
2[l00]  when  a  resolved  shear  stress  of  at  least  0.2  GPa  was 
applied.  Therefore,  the  core  with  full  Burgers  vector,  b  =  [100], 
is  lattice  trapped  and  is  essentially  unstable.  The  resolved  stress 
applied  to  assist  the  split  of  the  dislocation  into  partials  is  close 
to  the  resolved  stress  required  to  move  the  b  =  1  /2  [  100] 
partial.  This  perturbation  helps  the  system  to  reduce  its  energy 
into  a  configuration  compatible  with  the  minimum  of  the  y- 


surface  mentioned  above.  Dislocations  of  the  (010)[00l]  type 
did  not  split  and  relaxed  to  narrow  core  structures  with  b  = 
[001  ] .  The  calculated  widths  of  all  stable  configurations  are 
reported  in  Table  1. 


Table  1.  Ranking  of  Slip  Systems  Resulting  from  the 
Atomistic  Models 


slip 

plane 

type 

(GPa) 

e% 

(GPa) 

tsf 

(GPa) 

tpn 

(GPa) 

w  (A) 

(010) 

1/2[100]/ 

edge 

0.050 

0.012 

0.151 

0.201 

29.48 

(010) 

1/2[100]/ 

screw 

0.129 

0.034 

0.151 

0.280 

21.68 

(001) 

[0 10] /screw 

0.458 

0.044 

0.000 

0.458 

17.08 

(001) 

1/2[100]/ 

screw 

0.170 

0.012 

0.308 

0.478 

20.29 

(021) 

1/2[0T2]/ 

screw 

0.356 

0.039 

0.142 

0.498 

10.20 

(Oil) 

[l  00]  /screw 

0.518 

0.064 

0.000 

0.518 

14.11 

(021) 

[l  00] /screw 

0.524 

0.037 

0.000 

0.524 

8.85 

(021) 

1/2[100]/ 

edge 

0.296 

0.037 

0.280 

0.576 

16.95 

(Oil) 

[ 100] /edge 

0.583 

0.064 

0.000 

0.583 

20.25 

(001) 

1/2[100]/ 

edge 

0.287 

0.048 

0.308 

0.595 

36.89 

(021) 

[012] /screw 

0.635 

0.039 

0.000 

0.635 

8.33 

(001) 

1/2[010]/ 

edge 

0.240 

0.044 

0.406 

0.646 

30.07 

(021) 

1/2[100]/ 

screw 

0.382 

0.037 

0.280 

0.662 

14.94 

(010) 

[00 l] /edge 

0.782 

0.059 

0.000 

0.782 

12.77 

(010) 

[00 1] /screw 

0.831 

0.059 

0.000 

0.831 

14.77 

(011) 

[Oil]  /screw 

0.896 

0.075 

0.000 

0.896 

6.51 

3.1.2.  (001)  Plane.  Experimental  studies8  have  suggested  the 
(001)  plane  to  be  a  cleavage  plane  in  RDX.  Molecular  dynamics 
simulations  of  shocked  RDX  crystals5  reported  nucleation  of 
0.16[010]  partials  on  the  (001)  planes.  Ductile  behavior  was 
reported  for  both  the  (00l)[l00]  and  (00l)[010]  slip  systems 
in  ref  13  with  stable  stacking  faults  at  b  =  1/2  [  100]  and  b  =  1/ 
2[010].  The  [010]  and  [100]  directions  in  this  plane  are 
investigated  in  this  study. 

Edge  and  screw  dislocations  with  b  =  [100]  split  into  partials 
of  1/2[100]  type  upon  the  application  of  a  resolved  shear 
stresses  of  at  least  0.4  GPa.  This  lattice  trapping  effect  is  similar 
to  that  discussed  for  the  (010)  plane.  The  width  of  (00l)[l00] 
dislocations  are  reported  in  Table  1. 

Only  screw  type  full  dislocations  with  b  =  [010]  were  stable. 
Edge  dislocations  with  the  same  Burgers  vector  were  unstable 
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Figure  5.  Core  structure  of  the  (001)  [010]  edge  dislocation  with  b  =  1/2[010] .  (a)  Relative  displacement  of  COM  across  the  glide  plane,  (b) 
Molecular  conformational  state  of  the  core  in  (a).  Caaa  molecules  are  colored  in  red,  Caee  are  in  blue,  and  Caae  are  in  gray. 


x/l012  x/l100 


Figure  6.  Core  structure  of  b  =  1/2  [100]  dislocations  residing  on  the  (021)  plane 


.  (a)  Screw  dislocation,  (b)  Edge  dislocation. 


and  only  partials  with  b  =  1/2[010]  could  be  stabilized.  The 
core  structures  of  the  (00l)l/2[010]  edge  dislocations  are 
shown  in  Figure  5a.  For  these  partials,  rotation  and  conforma¬ 
tional  flips  (Caae  to  Caee)  of  the  molecules  across  the  glide 
plane  are  observed  and  an  out-of-plane  component  with  b  = 
0.1  [100]  develops. 

The  atomistic  core  structure  of  the  edge  partial  is  shown  in 
Figure  5b.  Conformation  changes  of  the  molecules  are 
observed  in  the  core  and  in  the  region  of  the  stacking  fault. 
The  molecules  colored  in  blue  are  in  the  Caee  state,  those  in 
red  are  in  the  Caaa  state,  and  all  others  are  in  the  Caae  ground 
state.  Most  molecules  in  the  stacking  fault  are  in  the  Caae 
conformation,  but  are  rotated  relative  to  the  molecules  in  the 
perfect  crystal.  Different  packing  is  observed  in  the  stacking 
fault  region.  It  should  be  also  noted  that  the  rotation  of  the 
molecules  leads  to  a  shift  of  the  COM  which  reflects  in  the  u(x ) 
profiles  shown  in  Figure  5a.  The  apparent  out-of-plane 
component,  «[ioo]0*)>  as  well  as  the  deviation  of  the  U[010](x) 
trace  from  the  Volterra  solution  in  the  stacking  fault  region 
(right  side  of  the  plot)  are  due  to  this  COM  shift  and  do  not 
represent  physical  dislocation  components  or  a  real  variation  of 
the  Burgers  vector  length. 

3. 1.3.  (021)  Plane.  The  (021)  plane  contains  the  slip  systems 
(02l)[l00]  and  (02l)[012].  Edge  and  screw  dislocations  with 
b  =  [100]  and  b  =  1/2[100]  were  found  to  be  stable.  This 
agrees  with  the  presence  of  a  minimum  in  the  y-surface 
reported  for  the  (02l)[l00]  system  at  the  1/2 [  100]  position.13 
Upon  application  of  a  resolved  shear  stress  of  0.55  GPa,  an 
edge  dislocation  with  the  full  Burgers  vector,  b  =  [100],  was 
found  to  lead  to  cleavage  starting  from  the  core,  which  indicates 


that  the  lattice  trapping  is  too  strong  in  this  direction  to  let  the 
core  relax  and  split  into  partials.  The  screw  dislocation  with  b  = 
[100]  does  not  lead  to  cleavage  but  is  also  incapable  of 
relaxation  by  splitting  into  partials.  The  partials  with  b  =  1/ 
2[l00]  are  stable  in  both  the  edge  and  screw  configurations, 
and  their  core  profiles  are  shown  in  Figure  6.  As  with  the  (001) 
[010]  partial  of  Figure  5,  apparent  out-of-plane  components  of 
b  =  0.l[012]  type  are  observed  upon  relaxation;  these  are  due 
to  COM  shifts  associated  with  molecular  conformation  changes 
in  the  core.  Screw  dislocations  with  b  =  [012]  and  b  =  1/ 
2  [012]  were  found  to  be  stable,  and  the  calculated  widths  are 
reported  in  Table  1. 

3.1.4.  (011)  Plane.  The  [100]  and  [Oil]  directions  are 
considered  in  this  plane.  The  possibility  of  cross-slip  in  the 
[100]  direction  was  postulated  in  previous  studies.12,13  The 
(01l)[01l]  slip  system  was  reported  to  exhibit  brittle  response. 
A  stable  stacking  fault  at  b  =  1/2[100]  was  reported  in  ref  13 
but  was  not  reproduced  in  our  atomistic  simulations.  The 
widths  of  stable  dislocations  with  full  Burgers  vector  on  the 
(01l)[l00]  and  (01l)[01l]  slip  systems  are  reported  in  Table 
1. 

3.2.  Ranking  of  Slip  Systems.  3.2.1.  Peierls  Stress  from 
Atomistic  Simulations.  The  stable  dislocations  are  ranked  in 
Table  1  in  the  order  of  increasing  tpn.  Table  1  also  includes  the 
critical  stress  applied  in  the  simulation  (ra)  which  puts  the 
dislocation  in  motion,  the  contribution  from  the  stacking  fault 
energy  (tSf)  to  the  Peach— Koehler  force,  and  the  width  ( w ) 
calculated  by  fitting  the  relaxed  core  structure,  as  described  in 
section  3.1. The  contribution  of  the  stacking  fault  to  the  critical 
stress  is  calculated  as 
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‘SF  ' 


?SF 

b 


(1) 


where  ySP  is  the  stable  stacking  fault  energy  and  b  is  the 
magnitude  of  the  Burgers  vector.  The  Peierls  stress  is  calculated 
as 


TPN  —  Ta  +  rSF  (2) 

indicating  the  fact  that  the  stacking  fault  aids  the  motion  of 
partials  in  all  models  considered.  The  error  in  ez,  is  the 
range  in  which  the  Peierls  stress  is  observed.  Specifically,  the 
dislocation  moves  at  the  stress  denoted  as  fa)  but  it  does  not 
move  at  ra  —  eT. 

The  Peierls  stresses  for  stable  dislocations  are  shown  in 
Figure  7  as  a  function  of  the  theoretical  strength.  The  bars  in 


Theoretical  Strength  (GPa) 

Figure  7.  Peierls  stress  values  evaluated  in  atomistic  simulations  versus 
the  theoretical  strength  in  the  respective  crystal  direction.  The  “e”  and 
“s”  labels  stand  for  edge  and  screw  dislocations,  respectively.  The  bars 
represent  the  range  in  which  the  Peierls  stress  is  observed.  The  limits 
defining  the  two  shaded  areas  are  selected  arbitrarily  to  correspond  to 
Gt/20  and  Gv/10,  where  Gv  is  the  Voigt  shear  modulus,  eq  4.  Gv/10  is 
close  to  the  theoretical  shear  strength  of  the  crystal. 


Figure  7  represent  eT.  The  theoretical  strength  is  calculated 
using  the  formula”1 


tts 


(3) 


Here  G  represents  the  shear  modulus  in  the  appropriate 
crystallographic  direction,  a  is  the  intermolecular  spacing  in  the 
b  direction,  and  d  is  the  interplanar  spacing  for  the  respective 
glide  plane.  One  observes  a  clear  demarcation  in  the 
distribution  of  the  Peierls  stress  values  close  to  350  MPa. 
This  value  is  close  to  Gv/20,  where  Gv  is  the  Voigt  average 
shear  modulus,  calculated  as22 

Gv  =  “(Qi  +  C22  +  C33  —  C23  —  C13  —  C12) 

+  -(C^  +  C55  +  CJ  ^ 


Here  C„  are  the  components  of  the  stiffness  matrix.  The 
calculated  value  of  Gv  at  0  K  using  the  moduli  reproduced  by 
this  potential  is  equal  to  7.004  GPa,  and  the  tpn  values  are 
discussed  here  in  terms  of  the  fraction  of  Gv. 

Dislocations  on  the  (010)[l00]  slip  system  have  the  lowest 
Peierls  stress  (rPN  <  Gv/20).  Thus  the  (010)  plane  is  found  to 
be  the  primary  slip  plane,  in  agreement  with  multiple 
experimental  observations.  For  slip  systems  of  the  type  (021) 


[012],  (02l)[l00],  (Oil) [100],  (00l)[010],  and  (00l)[l00], 
the  simulations  yield  Gv/20  <  tpn  <  Gv/10.  The  third  regime 
consists  of  slip  systems  of  the  type  (Oil)  [Oil]  and  (010)[00l], 
with  tpn  >  Gv/10.  Out  of  the  slip  systems  in  the  third  regime, 
(OlO) [OOl]  systems  can  be  considered  inactive,  since  the 
critical  stresses  calculated  are  close  to  the  theoretical  strength 
and  cleavage  is  likely  to  happen  before  plastic  deformation.  This 
result  is  not  in  agreement  with  previous  suggestions  made 
based  on  indirect  experimental  observations.  ’ 

It  is  noted  here  that  screw  dislocations  of  the  [100]  type  on 
the  (021),  (OOl),  and  (Oil)  planes  have  critical  stress  values  of 
0.524,  0.478,  and  0.518  GPa,  thus  indicating  the  possibility  of 
cross-slip  between  these  planes,  in  situations  where  the  resolved 
stresses  are  higher  than  0.524  GPa.  Cross-slip  has  been 
suggested  between  planes  containing  the  common  [100]  slip 
direction  in  other  studies.12'13,15  The  present  results  confirm 
that  this  is  possible  from  the  point  of  view  of  the  lattice 
resistance  to  slip.  The  unstable  stacking  fault  energies  for  these 
slip  systems  reported  in  the  literature13  are  also  within  a  narrow 
range  of  each  other.  The  critical  stress  for  the  [100]  screw 
dislocation  on  the  (010)  plane  is  0.28  GPa,  which  is  much 
lower  than  the  values  for  the  other  planes. 

3.2.2.  Peierls  Stress  Calculations  Using  the  Peierls- 
Nabarro  Model.  The  Peierls— Nabarro  (PN)  model  is  based 
on  an  interpretation  of  the  core  of  a  dislocation  as  a  continuous 
distribution  of  infinitesimal  dislocations  on  the  glide  plane.  The 
relevant  parameter  describing  the  core  shape  becomes  the 
density  of  these  infinitesimal  dislocations,  p  -  du/dx.  The 
model  requires  the  equilibrium  of  the  distribution,  which  is 
written  in  terms  of  the  force  resulting  from  the  elastic 
interaction  of  the  infinitesimal  dislocations,  the  lattice  rebound 
force,  and  the  applied  far-field  resolved  shear  stress.  All  PN 
models  assume  that  the  core  is  planar.  The  equilibrium 
condition  is  written  as 


G_  1  d (»(*')) 

2/ru_co  x  —  x'  dx' 


dx’  =  Fy(u(x))  +  t 


where 


(5) 


G' 


r\  sin2  0 

G[r^ 


and 


dy(ti) 

d  u 

(6) 


Here  G I  In  is  the  energy  prefactor,  9  is  the  angle  between 
vectors  b  and  £,  G  is  the  isotropic  shear  modulus,  v  is  the 
Poisson’s  ratio,  Fy  is  the  lattice  restoring  force  obtained  as  the 
derivative  of  the  y-surface  trace  in  the  direction  of  b,  and  r  is 
the  applied  resolved  shear  stress.  A  normalization  condition  is 
required,  which  is  written  as 


d(u(x')) 

dx' 


dx'  =  b 


(7) 


The  equilibrium  core  structure,  u(x),  is  obtained  by  solving  eq 
5  without  the  applied  shear  stress.  The  unknown  function  u(x) 
is  written  as  a  superposition  of  arctan  functions  as 


(8) 


The  Peierls  stress  can  be  obtained  by  increasing  the  value  of  r 
until  no  solutions  are  found  for  eq  5.  This  value  of  r  is  the 
critical  stress  required  to  move  the  dislocation  on  the  glide 
plane  (tpn)-  Note  that  this  version  of  the  PN  model  provides 
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only  an  approximation  of  the  dislocation  core  since  it  does  not 
include  nonlocal  effects  and  out-of-plane  displacements 
(shear— tension  coupling).  If  the  relative  motion  of  the  crystal 
planes  involves  dynamic  changes  of  the  conformation  state  of 
the  molecules,  the  PN  model  cannot  account  for  this  physics 
even  when  the  restoring  force  is  computed  from  the  y-surface 
evaluated  from  atomistic  models  (eq  6). 

y-surfaces  for  the  slip  systems  considered  in  this  study  have 
been  published  elsewhere.13  The  Peierls  stress  values  obtained 
from  the  PN  model  (t£n)  are  compared  with  the  values 
obtained  from  atomistic  simulations  in  Figure  8. 


Figure  8.  Comparison  between  Peierls  stress  values  obtained  from  the 
PN  model  (r(N)  and  atomistic  simulations  (rPN). 


The  (010)[l00]  system  is  found  to  be  the  most  easily 
activated  system  with  tpn  <  GJ 20  in  agreement  with  the 
predictions  obtained  from  atomistic  simulations.  Slip  systems  of 
the  types  (010)[00l]  and  (Oil) [Oil]  and  edge  dislocations 
with  b  =  1/2[010]  on  (001)  planes  were  found  to  have  t(n  > 
Gv/10.  The  values  obtained  from  the  PN  model  for  all  slip 
systems  are  within  20%  of  the  atomistic  results,  except  for  the 
(021) [0l2],  (Oil) [100],  and  (00l)[l00]  slip  systems  for 
which  the  error  is  about  25%.  This  provides  support  for  the  use 
of  the  PN  model  when  estimating  Peierls  stresses.  The 
computational  effort  associated  with  estimating  the  y-surface 
and  solving  the  PN  model  for  a  particular  slip  system  is  much 
smaller  than  that  required  to  determine  the  Peierls  stress  from 
full  atomistic  simulations,  since  very  large  atomistic  models  are 
required  in  the  latter  case. 

4.  DISCUSSION 

The  atomistic  results  can  be  analyzed  from  a  classical 
continuum  perspective.  Compatibility  conditions  require  that 
five  independent  slip  systems  are  required  in  a  crystal  to 
accommodate  an  arbitrary  applied  strain.  This  constraint  is 


relaxed  in  complex  ways  near  surfaces.  The  ability  of  the  RDX 
crystal  to  accommodate  an  arbitrary  plastic  strain  can  be 
investigated  using  the  method  proposed  by  Kocks23  and  Kelly 
and  Groves24  for  low  symmetry  crystal  structures.  This  method 
was  recently  employed  for  RDX  using  information  resulting 
indirectly  from  nanoindentation  experiments.12  The  strains 
produced  by  unit  slip  on  all  slip  systems  considered  in  this 
study  are  shown  in  Table  2. 

Assuming  that  all  eight  slip  systems  considered  are  activated, 
one  can  form  a  total  of  56  (SC5)  unique  combinations  of 
possible  groups  of  five  slip  systems.  Though  this  is  the  case, 
none  of  these  56  groups  possess  the  required  five  independent 
slip  systems.  At  the  most,  one  obtains  four  independent 
systems  and  this  points  to  the  fact  that  plasticity  in  RDX  is 
hindered  even  when  all  eight  slip  systems  investigated  are 
active. 

Additional  considerations  come  into  play  when  the  ranking 
of  slip  systems  presented  in  Table  1  is  taken  into  account. 
Recent  nanoindentation  experiments12  performed  on  RDX 
reported  first  pop-in  events  at  246,  368,  and  994  ^N  applied 
indentation  force  for  indents  in  the  (210),  (02l),  and  (001) 
faces,  respectively.  A  simple  analysis  based  on  Schmid  factors 
shows  that,  without  taking  into  account  the  Peierls  stress,  five, 
three,  and  two  slip  systems  could  be  activated  for  the  (210), 
(021),  and  (001)  cases,  respectively.  Let  us  consider  the  point 
located  at  0.48a,  with  a  being  the  contact  radius  of  the  indenter, 
below  the  surface  and  on  the  symmetry  axis.  In  an  isotropic 
material,  this  is  the  point  where  the  maximum  shear  stress  is 
reached  and  is  also  the  point  where  homogeneous  dislocation 
nucleation  is  usually  assumed  to  take  place.  We  compute  the 
resolved  shear  stress  in  all  slip  systems  studied  with  the  goal  of 
estimating  the  number  of  systems  in  which  activity  may  take 
place  at  the  critical  indentation  loads  observed  in  experiments. 
These  resolved  shear  stresses  are  plotted  against  the  lowest 
values  of  rPN  for  the  respective  slip  systems  in  Figure  9. 

For  indentation  on  the  (210)  face,  activity  is  enabled  on  the 
(010)[l00]  system,  for  indentation  on  the  (021)  plane  only 
screw  dislocations  on  (00l)[010]  can  move,  and  indentation 
on  the  (001)  face  leads  to  the  activation  of  (02l)[012]  at  the 
load  corresponding  to  the  first  pop-in  event.  This  indicates  that 
at  the  first  plastic  event  only  one  slip  system  is  activated  in  all 
situations.  This  is  meaningful  since  one  expects  that  the  first 
sign  of  plasticity  under  the  indenter  may  be  produced  by 
activity  on  one  slip  system  only.  However,  the  critical  loads 
predicted  using  our  data  would  underestimate  the  loads 
observed  in  experiments.  Specifically,  the  (010)[l00]  system 
is  activated  when  the  applied  force  normal  to  the  (210)  face  is 
12  ;/N,  the  (00l)[010]  system  would  operate  at  a  force  normal 
to  the  (021)  face  of  208  ^N,  and  the  (01l)[01l]  system 
becomes  active  at  a  critical  load  perpendicular  to  the  (001)  face 


Table  2.  Components  of  Strain  Resulting  from  Unit  Slip  on  Specific  Slip  Systems 


slip  system 

slip  plane  normal 

slip  vector 

£ll 

£22 

e33 

2^23 

2^31 

2^12 

(010)[100] 

[0  1  0] 

[10  0] 

0 

0 

0 

0 

0 

1 

(010)[001] 

[0  1  0] 

[0  0  1] 

0 

0 

0 

1 

0 

0 

(001)[100] 

[0  0  1] 

[10  0] 

0 

0 

0 

0 

1 

0 

(001)[010] 

[0  0  1] 

[0  1  0] 

0 

0 

0 

1 

0 

0 

(021 ) [100] 

[0  0.877  0.481] 

[10  0] 

0 

0 

0 

0 

0.481 

0.877 

(021 )  [0l2] 

[0  0.877  0.481] 

[0  -0.481  0.877] 

0 

-0.422 

0.422 

0.538 

0 

0 

(011)[100] 

[0  0.674  0.739] 

[10  0] 

0 

0 

0 

0 

0.739 

0.674 

(on)[oTi] 

[0  0.674  0.739] 

[0  -0.739  0.674] 

0 

-0.498 

0.498 

-0.092 

0 

0 

5332 
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Figure  9.  Activation  of  slip  systems  in  indentation  experiments 
perpendicular  to  the  (001),  (021),  and  (210)  crystal  planes,  “e”  and 
“s”  stand  for  edge  and  screw  dislocations,  respectively. 


of  295  ^/N.  The  present  analysis  predicts  the  proper  anisotropy 
of  the  indentation  response,  with  the  (210)  indents  having  the 
lowest  critical  load  and  the  (001)  indents  having  the  highest 
one.  The  underestimation  is  also  expected.  If  the  pop-in  event 
is  produced  by  dislocations  that  nucleate  homogeneously  under 
the  indenter,  the  load  required  for  nucleation  should  be  much 
larger  than  that  required  for  motion  of  an  existing  dislocation 
since  the  respective  energetic  barriers  are  quite  different. 
Likewise,  if  the  pop-in  event  is  due  to  preexisting  dislocations, 
the  critical  load  should  not  only  put  these  dislocations  in 
motion,  but  should  also  be  large  enough  to  render  unstable  the 
preexisting  dislocation  structures.  Therefore,  in  both  situations, 
the  first  pop-in  loads  predicted  based  on  the  Peierls  stress 
values  are  expected  to  be  smaller  than  those  observed 
experimentally.  Since  the  mechanism  leading  to  the  first  pop- 
in  is  not  entirely  known,  the  observation  made  here  that  the 
ranking  of  the  activated  slip  systems  in  terms  of  Peierls  stresses 
matches  the  experimental  anisotropy  in  indentation  should  not 
be  overemphasized.  We  interpret  it  as  a  suggestion  that,  since 
Peierls  stresses  are  very  large  in  these  crystals,  the  critical  stress 
required  to  put  dislocations  in  motion  should  be  close  to  the 
critical  stress  of  the  mechanism  leading  to  criticality  in 
indentation. 

Let  us  consider  the  implications  of  the  results  presented  in 
this  article  in  the  context  of  anisotropic  sensitivity  of  RDX 
shocked  perpendicular  to  the  (ill),  (210),  and  (100)  planes. 
Calculation  of  Schmid  factors  shows7  that  (ill)  shocks  could 
activate  seven  slip  systems,  (210)  shocks  provide  nonzero 
resolved  shear  stress  on  five  slip  systems,  and  (100)  shocks  do 
not  activate  any  slip  system,  even  when  the  (001)  systems  are 
taken  into  consideration,  (ill)  shocks  activate  systems  (001) 
[100],  (00l)[010],  and  (010)[00l],  which  are  not  available  for 
(210)  shocks. 

Interesting  observations  can  be  made  by  considering  specific 
shock  strengths.Pj^j  values  of  magnitude  7,  2.14,  and  1.24  GPa 
are  considered.  These  values  are  chosen  on  the  basis  of 
previous  experimental  studies4’5,7  and  represent  critical  thresh¬ 
olds  for  RDX.  Pj^  =  7  GPa  is  the  minimum  experimental  peak 
stress  reported  for  decomposition  in  RDX.7  Low  stress  flyer- 
plate  impact  experiments  reported  an  abrupt  change  in  the 
response  of  RDX  for  (ill)4  and  (02l)6  shocks  for  1.24  <  Prh 
<  2.14  GPa.  Values  of  the  resolved  shear  stress  (trss)  f°r  all  slip 
systems  are  shown  in  Figure  10  against  the  lowest  values  of  rPN 
in  the  respective  system,  for  (ill)  shocks  with  Pj^  =  7,  2.14, 
and  1.24  GPa  and  for  a  (210)  shock  withF^  =  7  GPa. 
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For  a  Phh  of  7  GPa  it  results  that  all  seven  slip  systems  with 
nonvanishing  Schmid  factors  have  Trsj  >  tpn  for  (ill)  shocks. 
(210)  shocks  activate  five  slip  systems.  However,  the  sensitivity 
is  higher  for  (210)  shocks.7  On  the  basis  of  this  observation  a 
lower  number  of  slip  systems  results  in  hindered  plasticity  and 
the  inability  of  the  crystal  to  relax  the  applied  stress.  Thus  the 
present  study  confirms  the  validity  of  the  steric  hindrance 
model  pertaining  to  the  decomposition  of  RDX.  A  recent 
computational  study  of  compressive  dynamic  loading  using  a 
reactive  force  field  also  supports  the  validity  of  the  steric 
hindrance  model.25  It  should  be  noted  that  the  RDX  crystals 
undergo  an  a— y  phase  transition  at  pressures  above  4.5  GPa26 
for  shock  and  between  2  and  4  GPa27'28  under  quasi-static 
loading  conditions.  In  addition  to  the  phase  transition,  the  high 
strain  rates  and  temperatures  involved  in  a  shock  might  result 
in  more  complicated  deformation  mechanisms  compared  to 
what  can  be  inferred  from  quasi-static  simulations. 

For  Prjj  =  2.14  GPa  in  a  (ill)  shock,  at  least  five  slip 
systems  (with  m  ^  0)  remain  active.  If  one  considers  ^rh  - 
1.24  GPa  and  the  (ill)  shock,  only  the  (010)[l00]  system  is 
clearly  activated.  Thus  the  observed  change  in  response  of 
(ill)  oriented  crystals  can  be  attributed  to  the  activation  of  a 
higher  number  of  slip  systems  as  the  shock  pressure  increases. 
Also,  for  (ill)  shocks,  the  cross-slip  of  [  1 00]  screw 
dislocations  will  be  inhibited  at  lower  shock  pressure  and  this 
might  also  play  a  role  in  controlling  the  deformation  behavior. 

In  the  case  of  (021)  shocks  (not  shown),  Pj^j  =  2.14  GPa 
activates  the  (00l)[010]  and  (010) [001]  slip  systems.  A  shock 
of  Ppjj  =  1.24  GPa  will  activate  only  the  (00l)[010]  system. 
Thus  the  above  conclusion  that  anisotropic  response  can  be 
correlated  with  differential  activation  of  slip  systems  remains 
unchanged  for  the  (021)  shocks. 

Shocks  normal  to  the  ( 100)  plane  activate  none  of  the  slip 
systems  considered  in  this  work.  It  is  possible  that  other 
deformation  mechanisms4,29  and/or  other  slip  systems  become 
active  in  order  to  accommodate  the  imposed  deformation. 

5.  CONCLUSIONS 

Atomistic  simulations  in  conjunction  with  the  Peierls— Nabarro 
model  were  used  to  study  the  mobility  of  dislocations  in  RDX. 
This  study  indicates  that  the  (010)[l00]  slip  system  has  the 
least  resistance  to  dislocation  motion.  Screw  dislocations  of  the 
[  100]  type  on  (021),  (001),  and  (Oil)  planes  were  found  to 
have  Peierls  stress  values  within  10%  of  each  other,  supporting 
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the  conjecture  that  cross-slip  of  screw  segments  with  [100] 
orientation  is  possible.  Critical  stresses  were  found  to  be 
between  Gv/20  and  Gv/10  for  the  majority  of  the  slip  systems 
investigated.  However,  the  number  of  independent  slip  systems 
is  always  lower  than  the  minimum  of  five  required  for  the 
accommodation  of  a  general  strain  state;  hence  plasticity  in 
RDX  is  hindered.  The  anisotropic  sensitivity  of  RDX  crystals 
shocked  in  (ill)  and  (210)  directions  can  be  explained  in 
terms  of  differential  mobility  on  available  slip  systems,  and  this 
supports  the  validity  of  the  steric  hindrance  model  for  this 
material.  An  explanation  for  the  anomalous  response  at  low 
shock  pressures  is  also  provided  in  terms  of  the  activation  of 
available  slip  systems. 
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